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ABSTRACT 

We present fully sampled 12 C0, 13 C0 J=2-l, 3-2 maps of the inner 
~ 1' x 1' region of NGC 1068. We combine these measurements with an 
existing interferometric map of 12 CO J =1-0 that includes single dish data and 
thus contains all the flux present. This allows a reliable estimate of the 12 CO 
(J=3-2)/(J=l-0) ratio at the highest angular resolution currently possible and 
the use of this sensitive line ratio to probe the physical condition of the molecular 
gas. We also present two measurements of the faint C 18 J=2-l emission 
which confirm earlier measurements of a high C 18 0/ 13 CO intensity ratio in 
this galaxy. The ratios of the 12 CO, 13 CO isotopes can only be reproduced for 
small/moderate optical depths of 12 CO J=l-0 (r ~ 1 — 2) which is incompatible 
with the high C 18 0/ 13 CO ratios observed. A simple two-phase model for the 
gas can account for all the observed line ratios if the C 18 emission and part 
of the 13 CO emission arise in a dense spatially concentrated component, where 
C 18 J=l-0 has optical depths of r ^ 1. The 12 CO emission originates from 
a warmer, diffuse gas phase with r ~ 1 — 2 for J=l-0. The dense gas phase 
contains the bulk of the molecular gas mass while the diffuse phase may not 
be virialized leading to an overestimate of molecular gas mass when deduced 
from the luminosity of the 12 CO J=l-0 line and a standard galactic conversion 
factor. This suggests that, since type 2 Seyferts harbor a central starburst more 
often than type 1, the higher average 12 CO J=l-0 luminosity of type 2 hinted 
by earlier studies may simply reflect a difference in molecular gas excitation 
rather than in gas mass. 

Subject headings: galaxies: individual (NGC 1068) — galaxies: Seyfert — galaxies: 
starburst — ISM: molecules 
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1. Introduction 

Intense central star formation and possibly an Active Galactic Nucleus (AGN) are 
thought to be the energy sources responsible for the observed large IR luminosities in bright 
IRAS galaxies (e.g., [Telesco 1988fc |Sanders, Scoville fc Soifer 1991|) where large amounts 



of molecular gas are found ( panders, Scoville fc Soifer 199 1|) . An important problem in 
understanding the energetic phenomena occurring in the nuclear regions of such galaxies is 
determining the role of the molecular gas in fueling the AGN and/or a central starburst. 

Many of the very IR-luminous (L FIR > 1O 11 L ) galaxies are found to be strongly 
interacting or merging systems QSanders et al. 1988fc [Lawrence" et al. 1989|) . Strong 



interactions between two galaxies, at least one of which is gas rich, may foster the rapid 
accumulation of large amounts of molecular gas in the nuclear regions (|Barnes fc Hernquist 



1991|) where they can "ignite" a starburst and "fuel" an AGN. For Seyfert galaxies with 
more moderate FIR luminosities (Lpm ^ 10 11 L Q ) the picture is less clear. Earlier work by 
Heckman et al. (1989) suggested that the hosts of type 2 Seyfert nuclei contain, on average, 
more molecular gas mass and have a larger FIR luminosity and hence higher star forming 
activity than type 1. More recent studies by Maiolino et al. (1997) of a larger and better 
defined sample of Seyferts and field spirals confirms the higher FIR luminosities of Seyfert 
2's but casts some doubts on the purported difference in molecular gas mass between these 
galaxies and Seyfert l's, or field spirals. In the same study preliminary evidence points 
towards bars, interactions and distorted morphologies as being present more frequently in 
the hosts of type 2 than in type 1 Seyferts. In this case, efficient molecular gas transport 
towards the nucleus will be present more often in former rather than the latter. Hence, 
provided that this mechanism remains effective over scales ranging from L ~ 1 kpc down to 
a few parsecs, it can naturally explain why the onset of a central starburst correlates with 
a higher probability of obscuring the Broad Line Region of the AGN. 
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Nevertheless such a scenario still cannot explain a larger total molecular gas reservoir 
in Seyfert 2 with respect to Seyfert 1 galaxies. However, since type 2 are more likely to 
harbor a central (~ 1 — 5 kpc) starburst than Seyfert l's (Maiolino 1995), there exists 
the possibility that any difference in their average 12 CO J=l-0 luminosities is caused by 
differences in the average excitation of the molecular gas rather than the total gas mass. 

NGC 1068 is the closest and best studied example of a luminous Seyfert 2 galaxy 
which also harbors an intense starburst in its central ~2.5 kpc. It has a large bolometric 
luminosity of ~ 3 x 10 11 L roughly equally divided between the unresolved AGN and the 
starburst disk ( [l'elesco et al. 1984] ) . The distribution of molecular gas is characterized by 
a nuclear concentration within the central 100 pc (Planesas, Scoville & Myers 1991; Heifer 
& Blitz 1995) where mostly dense gas (n(H 2 ) ~ 10 5 cm~ 3 ) resides (Jackson et al. 1993; 
Tacconi et al. 1994; Heifer & Blitz 1995) and a more extended molecular gas distribution 
that coincides with the starburst region and consists of two inner spiral arms and a bar 
(Heifer & Blitz 1995). The nuclear component is thought to contribute to the built-up of 
an AGN-obscuring molecular torus (Planesas et al. 1991; Cameron et al. 1994; Tacconi et 
al. 1994). In the present work we use this galaxy as a "laboratory" to study the excitation 
of the bulk of the molecular gas in the intense starburst environment of the host galaxy 
of a Seyfert 2 nucleus. To do so we use our fully sampled maps of 12 CO, 13 CO J=3-2, 
2-1 acquired with the James Clerk Maxwell Telescope (JCMT)Q in combination with an 
existing interferometric 12 CO J=l-0 map (Heifer fc Blitz 1995| ). 

Resolution-matched observations of different rotational lines of 12 CO and especially 
13 CO provide an excellent probe of the physical conditions in the bulk of the molecular 

2 The JCMT is operated by the Joint Astronomy Center in Hilo, Hawaii on behalf of 
the parent organizations PPARC in the United Kingdom, the National Research Council of 
Canada and the The Netherlands Organization for Scientific Research. 
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gas. A ratio like (3-2)/(l-0) can be particularly sensitive to gas excitation conditions since 
the J=3 level lies 32 K above ground versus 5.5 K for J=l while the J=3-2 transition 
has a significantly higher critical density (n cr = 4 x 10 4 cm -3 ) than the J=l-0 transition 
(n cr = 2 x 10 3 cm -3 ) in the optically thin regime. Moreover the maps of the rare 13 CO 
isotope allow the study of the molecular gas using less opaque lines than the ones of 12 CO 
and therefore sensitive to the average conditions of a larger fraction of molecular material. 

Throughout this study we assume a distance to NGC 1068 of D=14 Mpc (Sandage & 
Tammann 1981, H Q = 75 km s _1 ) for which 1" corresponds to 68 pc. 

2. OBSERVATIONS 

We used the receivers A2, B3i and the newly commissioned receiver B3 on JCMT 
in several observing runs in order to produce fully sampled maps (Nyquist sampling) of 
the emission lines 12 CO, 13 CO J=3-2, 2-1 of the inner ~ 1' x 1' region of NGC 1068. We 
also observed the very weak C 18 J=2-l transition in two points, the details of all our 
observations are summarized in Table 1. 

EDITOR: PLACE TABLE 1 HERE 

The map center was located at a(B1950) = 02 h 40 m 07 s .06, 5(B1950) = -00° 13' 31.45" 
and the receiver was tuned to a velocity of vlsr = 1125 km s^ 1 . For more efficient mapping 
of the CO emission, the area mapped was along a rotated grid with PA = 40°. The grid-cell 
is an equilateral triangle with a size of A9 S ~ #hpbw/2, where 6*hpbw is the HPBW of the 
gaussian beam. This grid provides the optimum coverage of the CO-emitting region of 
NGC 1068 at a Nyquist sampling and is uniquely determined by the demand of a) a uniform 
sampling pattern and, b) every sampled point having the maximum number (N=6) of the 
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closest neighboring points at a distance of A6 S . We performed all our observations using 
beam switching with a beam throw of 150" in azimuth at the recommended rate of 1 Hz 
and only for the very weak C 18 J=2 1 the rate was 2 Hz in order to achieve a stable 
baseline. The pointing and focus were monitored frequently by observing bright quasars, 
planets and OMC1 (for 13 CO J=3-2). 

We converted the temperature scale of the JCMT spectra to the T mb scale by 
using the relation T mb = T* A /r] inh . We performed measurements of r] mh for all the receivers 
by observing Saturn (A2, B3i) and Jupiter (B3). We found excellent agreement (within 
~ 10%) with the values reported in the JCMT manual (Matthews 1996). The adopted 
beam efficiencies are 77 mb (A2) = 0.69, r] mb (B3i) = 0.55 and ?? mb (B3) = 0.62. The error 
associated with the line intensity measurements was estimated using the relation 



The first two terms express the thermal rms error and the stochastic uncertainty of 
the calibration of the line intensities. The third term is the systematic uncertainty of the 
assumed telescope efficiency factors. If T is the main beam brightness temperature averaged 
over N ch channels and N bas is the total number of channels defining a symmetric baseline 
around the spectral line, then 



where ST^ is the thermal noise per channel. 

The calibration uncertainty factor was estimated from the dispersion of the observed 
intensities of many spectral line standards with high S/N, and it was found to be 
~ 0.10 — 0.15 for all the observed frequencies. The third factor is ~ 0.10 for all frequencies 




(1) 
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(Friberg, Sandell, private communication). Finally we performed an additional check on 
the overall uncertainty of the line intensities by frequently monitoring the 12 CO lines in 
some standard points of our map. For the weak 13 CO lines we monitored the spectrum of 
OMC1 at Vlsr = 1125 km s _1 where various bright lines are still present. In all cases we 
found excellent agreement of the line profiles and a measured dispersion of line intensities 
consistent with the one estimated from the first two terms of Equation 1. 

3. DATA REDUCTION 

The JCMT data were reduced using the SPECX reduction package. We smoothed all 
our spectra to a common velocity resolution of AV c han ~ 8 km s _1 and removed linear 
baselines. In the case of the 12 CO J=3-2 transition we also removed some bad channels 
that were present at the high velocity end of the band. 

The grid maps were interpolated by convolving them with a gaussian function. For the 
J=3-2 maps the FWHM of the interpolating function used was Q\ = 8". The interpolated 
map has effective resolution of [6* i 2 + 6 |2 I p BW (3 — 2)] 1//2 = 16". For the J=2-l grid maps we used 
6i = 11" and the effective resolution of the interpolated map is [9f + #hpbw(2 — -Q] 1 ^ 2 ~ 24". 
The maps were further analyzed with the AIPS reduction package. The task HGEOM was 
employed to re-grid the JCMT maps in order to re-orient them in an (RA, DEC) grid. A 
check of the consistency of the interpolation and re-gridding was performed by comparing 
the area-averaged T m b in individual channel maps between the original grid maps and 
the final ones used in our analysis. We found good agreement within the thermal rms 
uncertainties. The maps of the velocity-averaged main beam brightness (T m b)Av for 12 CO, 
13 CO are shown in Figures 1 and 2 for the ,1=2-1 and 3-2 transitions respectively. 



EDITOR: PLACE FIGURE | HERE. 
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EDITOR: PLACE FIGURE | HERE. 

The 12 CO J=l-0 map used in our analysis was obtained by Heifer & Blitz (1995) with 
the Berkeley-Illinois-Maryland Association (BIMA) interferometer and the NRAO 12m 
telescope. This map is the most suitable one for comparison with our JCMT maps since 
its combination of interferometer and single dish measurements recovers all the CO flux. 
The interferometer map has been corrected for primary beam attenuation by assuming 
a gaussian profile with HPBW=100". We convolved the 12 CO J=l-0 channel maps to a 
spatial resolution of 16" and converted the original brightness units (mJy/beam) to the T m b 
temperature scale in order to compare them to the 12 CO J=3-2 map. The channel maps 
and the associated line ratio T32 = (3 — 2)/(l — 0) are shown in Figure 3. 

EDITOR: PLACE FIGURE g HERE. 

The close correspondence of the 12 CO J=3-2 and the 12 CO J=l-0 emission in every 
channel demonstrates that no serious pointing offsets are present between the two maps. 
The estimated channel-to-channel rms dispersion of the Y32 ratio is of the order of 10%. 
This includes only the thermal rms error since the other two sources of error for the JCMT 
maps (Equation 1) and the systematic calibration uncertainty in the flux conversion factor 
of the BIMA 12 CO J=l-0 map are constant across the passband. The largest single source 
of error is the ~ 30% flux calibration uncertainty of the 12 CO J=l-0 BIMA map. This 
results to a total uncertainty of 35% for the r 32 ratio. 

The C 18 J=2 1 spectra correspond to a spatial resolution of #hpbw — 22". In order 
to compare them to the 13 CO spectra towards the same locations we used an interpolated 
13 CO J=2-l map with Oq = 8", which gives an effective resolution of Q\ = 23". The spectra, 
smoothed to a common frequency resolution of 25 MHz (34 km s _1 ) are shown in Figure 4. 

EDITOR: PLACE FIGURE | HERE. 
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4. RESULTS 

The 13 C0, 12 C0 ,1=2-1 and 3-2 maps (Figures 1, 2) demonstrate that the relative 
intensities of these two isotopes change significantly across the emitting region. There is 
a significant increase of the intensity of 13 CO relative to 12 CO towards the location of 
the giant molecular gas associations (GMAs) in the southern part of the central starburst 
region of NGC 1068. These molecular gas concentrations show up in high resolution maps 
(Planesas, Scoville & Myers 1991; Heifer & Blitz 1995) as high brightness regions with large 
inferred H 2 surface densities. Especially large variations of 12 CO/ 13 CO are observed for 
the J=l-0 transition (Heifer & Blitz 1995; Papadopoulos, Seaquist & Scoville 1996). Large 
variations of this ratio have been recently reported by Aalto et al. (1997) in the extreme 
starburst /merger system Arp 299 where high resolution interferometric measurements 
indicate changes by a factor of ~ 6. 

The r 32 ratio also shows a large range of values in the various regions in the inner 
3 kpc of NGC 1068 (Fig. 3). The lowest values r 32 ~ 0.3 — 0.4 occur towards the eastern 
part of the emission while the highest ones r 32 ~ 0.6 — 0.8, towards the western and 
southwestern part where the bright GMAs are located. Moreover, the ratio r 2 i varies in a 
similar fashion with a range of r 2 i ~ 0.6 — 0.7 and r 2 i ~ 0.7 — 0.8 for the eastern and the 
western/southwestern parts of the emission respectively. 

Line ratios as low as r 32 ~ 0.3 — 0.4 characterize the bulk of the relatively cold gas in 
the Giant Molecular Clouds (GMCs) in the disk of the Galaxy flSanders et al. 1993| ). The 
same is true for r 21 ^0.7, usually measured in the outer disk of the Galaxy ( |Hasegawa 



19971) , or the disks of other galaxies (e.g., |Eckart et al. 1990| ; Eckart et al. 1991). An LTE 



approximation of such low ratios yields a Tkin ^ 6 K, comparable to the minimum value 
of 7 K permitted by cosmic ray heating (Goldsmith & Langer 1978). Then a maximum 
brightness of 12 CO J=l-0 would be Tb ~ 3 K, much smaller what is usually measured in 
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NGC 1068 ( |Plancsas, Scovillc fc Myers 1991|) , which implies sub-thermal excitation of the 
J=3-2, 2-1 transitions towards the regions with such low r 32 and r 2 i ratios. 

On the other hand ratios as high as r 32 ~ 0.6 — 1 are measured towards central regions 
of starbursts QDevereux et al. 1994| |Wild at al. 1992| ) and close to the cores of star forming 



GMCs like W3 in the Milky Way (e.g., [Phillips et al. 198 1| ). A ratio of r 21 ~ 0.80 - 1.0 



is found mainly towards the centers of nearby spirals ( [Braine fc Combes 1992|) , centers of 



starbursts (e.g., [Eckart et al. 1990| ; |Eckart et al. 1991| ; |Wild at al. ~1M% ), close to HII 



regions and star forming cores in the Galaxy ([Hasegawa 1997|) . In NGC 1068 the largest 



values (r 3 2 = 0.8 and r 2 i = 1) are observed towards the location of the most massive GMAs 
identified in the interferometer maps. Under LTE conditions these ratios imply average 
kinetic temperatures of Tk in ^ 25 K. In contrast, the bulk of the molecular gas in the plane 
of the Milky Way is significantly colder at Tk in « 10 K ( [Sanders et al. 1993] ), undoubtly the 
result of a more quiescent environment than the central starburst region of NGC 1068. 



4.1. The isotope line ratios 

As reported previously the value of 12 CO/ 13 CO ratio changes across the central kpc-size 
starburst region of NGC 1068 with the largest variations measured for the J=l-0 transition. 
This shows most prominently in high resolution interferometer maps (Papadopoulos et al. 
1996) but even in the area-averaged spectra shown in Figure 5. 

EDITOR: PLACE FIGURE | HERE. 

These spectra demonstrate that the ratio Rio = 12 CO/ 13 CO (J=l-0) changes by a 
factor of ~ 2 along the line profile while variations of R 2 i (J=2-l) and R 32 (J=3-2) follow 
in a similar fashion but are significantly less prominent. It is important to note that the 
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variation of Rio across the line profile is not influenced by the systematic uncertainties of 
the line intensities, since they can only introduce an overall shift of the T mb scale. The line 
shapes are subject only to the thermal rms error across the bandpass which is of the order 
of ~ 10% for the J=l-0, and ~ 15% for the J=2-l and 3-2 transitions. 

The reason for using the OVRO 13 CO J=l-0 spectrum (Papadopoulos et al. 1996), 
which does not include single dish data, instead of the one from BIMA/NRAO 12m (Heifer 
& Blitz 1995) stems from the fact the later has a significantly lower S/N, a result of poor 
weather during the observations (Heifer, private communication). Nevertheless both data 
sets exhibit the same basic characteristics. Moreover, the largest variation of Rio is found 
among the low and the high velocity end of the J=l-0 line profile, where the addition of 
the single flux data is not expected to alter the line intensity (Heifer & Blitz 1995). 

In the inner kiloparsec-size region of NGC 1068 we find Rj+i j ^ 10 (J=0, 1, 2) 
throughout, irrespective of location or averaging scale. For the entire range of the 
[ 12 CO/ 13 CO] abundance values (e.g., Langer & Penzias 1993; Henkel & Mauersberger 
1993) such high isotopic ratios correspond to low 12 CO optical depths. Assuming LTE 
and [ 12 CO/ 13 CO] = 40, a value of Ri ;> 10 corresponds to a 13 CO J=l-0 optical depth 
of Tio^ ^0.1, and hence rather high = r2i( 13 CO) and r^ = r32( 13 CO) ratios (or 
equivalently, high R 2 i, R32 ratios) with respect to the ones observed. This has been used 
to support the notion of two-phase gas in starburst centers (Aalto et al. 1995), yet there 
is still a wide range of non-LTE conditions that must be explored before such a conclusion 
can be reached. We will see that the best evidence for such a state of the molecular gas 
comes from the C 18 J=l-0 (Papadopoulos et al. 1996) and J=2-l data (this work). 
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4.2. The high C 18 0/ 13 CO ratio 

The ratio R^ 8) = C 18 0/ 13 CO J=2-l obtained at the two positions is Rg 8) = 0.26 ±0.07 
and R^ 8 "* = 0.20 ± 0.06 (Figure 4). These values, while somewhat lower, are comparable to 
the average ratio Rio 8 "* = 0.3 ±0.1 found from high resolution OVRO maps (Papadopoulos 
et al. 1996) for the J=l-0 transition. Moreover, these maps also show that the C 18 
emission does not smoothly trace 13 CO in the various GMAs and varies drastically from one 

(18) 

location to the next with values as high as Ri ~ 0.65. This is significantly higher than 
the abundance ratio [C 18 0/ 13 CO] = 0.12 — 0.14 inferred for the Milky Way from various 
studies (e.g., |Wannier 1980| ; |Hcnkel fc Mauersberger 1993| ; |hanger fc Penzias 1993j ). Ratios 



R-iQ 8 '* ^ 0-2 are also found in several central starbursts of galaxies (e.g., [Sage, Mauersberger 



fe Henkel 199g ; |Casoli, Durpaz, fc Combes 1992| ; |Henkel fc Mauersberger T993|) 



While it is possible that an enhancement of the [C 18 0/ 13 CO] abundance can occur in 
starburst environments (Henkel & Mauersberger 1993) and hence play a role in the observed 

(18) (18) 

high values of R 10 and R^ , it is very unlikely that it can be solely responsible for the 

high Riq 8 ' ) values observed in individual GMAs (Papadopoulos et al. 1996). Furthermore, 

(is) 

the large spatial variations of the K 10 ratio over scales of ~ 400 pc are difficult to attribute 
solely to a varying [C 18 0/ 13 CO] abundance because the mixing of interstellar gas is expected 
to be effective enough to homogenize the isotope ratios over similar or larger scales, even 
during an ongoing burst of star formation (Henkel & Mauersberger 1993). 

Thus, it seems more likely that the optical depth of the rare isotopes is not negligible. 
For RSo 8) ;> 0.20, the optical depth of 13 CO J=l-0 is r x ( J 3) £ 1.5 (LTE assumed), which 
directly contradicts the general conclusion reached earlier from the Rio ratio which yields 
Tio ^ 0.1. This constitutes the strongest evidence that 12 CO, 13 CO and C 18 cannot be 
all tracing the same phase of the molecular gas. 
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5. DISCUSSION 

Line ratios of CO and its isotopes as inputs to a radiative transfer model that solves the 
rate equations for the various rotational levels are useful probes of the physical conditions 
of the molecular gas. Obviously the larger the number of the input ratios the better the 
constraints on the gas properties. The various spectral lines for NGC 1068 have been 
observed at very different angular resolutions with the lowest one (6*hpbw ~ 24" for J=2-l) 
being comparable to the source size. Hence, using the spatially-averaged spectra over the 
entire CO emitting area (~ 1' x 1') allows the estimate of all the line ratios for that area 
and thus offers the maximum number of inputs that can be used to constrain the physical 
properties of the molecular gas. These line ratios are tabulated in Table 2. 

EDITOR: PLACE TABLE 2 HERE 

For the modeling of the molecular gas properties we used a Large Velocity Gradient 
(LVG) code. We searched an extensive grid of LVG models with a parameter range of 
T kin = 10 - 120 K (AT kin = 2 K), n(H 2 ) = (0.1 - 10 4 ) x 10 3 cm- 3 (Alogn(H 2 ) = 0.5) and, 
A = X/(dV/dr) = (0.1 - 10 2 ) x lO" 6 (km s" 1 pc" 1 )" 1 (AlogA = 0.5), where X = [ 12 CO/H 2 ] 
and dV/dr is the velocity gradient for the "average" molecular cloud. We assumed an 
abundance ratio of [ 12 CO/ 13 CO] = 40 since studies of the Milky Way (Langer & Penzias 
1993) and central regions of starburst galaxies (Henkel at al. 1993; Henkel & Mauersberger 
1993) show this value to be a suitable average of the abundances observed towards galactic 
centers. The abundance of [C 18 0/ 13 CO] = 0.15 is derived from [ 12 CO/ 13 CO] = 40 and 
[ 12 CO/C 18 0] = 250, which is the highest such value derived for molecular clouds in the 
Milky Way (Wannier 1980), and was measured towards its center. We examine the entire 
grid of the LVG models for the values of r 2 i, r3 2 that correspond to 0%, —30% and +30% 
offset in the BIMA 12 CO J=l-0 intensity. The best fit is found for the highest r 2 i, r 32 ratios 
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which correspond to the 30% reduction in the flux density of the 12 CO J=l-0 map. 

Apart from giving the best fit, such an offset of the 12 CO J=l-0 intensity observed 
by BIMA seem to be corroborated by single dish measurements reported by Maiolino 
et al. (1997) and Kaneko et al. (1989). The first used the NRAO 12m to obtain a 
velocity-integrated main beam brightness for the inner 55" x 55" of NGC 1068 which is 
I = (90 ± 13)K km s -1 , while the latter used the NRO 45m telescope to map a similar 
area and found I = (83 ± 17)K km s _1 . These values were obtained from their data after 
converting the T^ (NRAO 12m) and T^ (NRO 45m) temperature scales to T m b- The 
integrated flux found over the same region from the BIMA map is I = (131 ± 40)K kms^ 1 . 
Clearly only a large offset of —30% can bring the BIMA value in agreement with the other 
two. Moreover, adopting such an offset yields r 2 i = 0.97, which in much better agreement 
with the value r 2 i = 1.1 reported by Braine & Combes (1992) over the same angular size. 

The conditions corresponding to the best fit are T kin = 20 K, n(H 2 ) = 10 4 cm" 3 and 
A = 10~ 6 (km s _1 pc -1 ) -1 . It is interesting to note that Thronson et al (1987) using 
mm/sub-mm continuum measurements deduce a dust temperature of T=20 K for the inner 
55" of this galaxy, which is similar to the gas temperatures we find from the LVG fits. 
This is expected for well mixed gas and dust where the 12 CO, 13 CO lines trace the same 
"average" conditions as the mm/sub-mm emission originating from dust. 

As expected from our previous discussion, the best LVG model cannot reproduce the 
observed high C 18 0/ 13 CO line ratios, since the deduced optical depth of the 12 CO J=l-0 
transition is r 10 ~ 1 — 2. In this case 13 CO, C 18 are very optically thin and the ratio 
of their intensities is equal to their assumed relative abundance. The moderate optical 
depth of the 12 CO J=l-0 transition was suggested by Aalto et al. (1995) as being one of 
the characteristics of a second, more diffuse and warm gas phase that surrounds a dense 
and more spatially concentrated one in the nuclear regions of starbursts (see also Wall et 



- 15 - 



al. 1993 for earlier work on spirals). This would naturally explain why both the average 
brightness of the 12 CO J =1-0 transition as well as the R 10 line ratio are more sensitive 
functions of velocity (and therefore of location within the galaxy) than the J=2 1, 3-2 
transitions (Figure 5) and the R21, R32 ratios respectively 

For T kin > 20 K it is 2.7r 10 < r 2 i < 4r 10 and 3r 10 < r 32 < 9t 10 (LTE), hence a moderate 
T10 ~ 1 — 2 corresponds to r J+ i j ~ 3 — 6 (J=l,2) or higher. This range of r w values 
allows the 12 CO J=l— transition to be significantly more sensitive than the J=2-l and 
J=3-2 to variations of gas column density and excitation conditions. Higher values t 10 (> 2) 
would render all three transitions optically thick and thermalized even at moderate/low 
(n(H 2 ) ^ 10 3 cm -3 ) gas densities because of radiative trapping, and hence they would be 
equally "insensitive" to the conditions of the molecular gas except the kinetic temperature. 
On the other hand, significantly lower values r 10 (^ 0.1) yield small/moderate optical depths 
for all three transitions and not just J=l-0. Hence they will all be sensitive to the gas 
excitation conditions and their intensities and Rj+i j line ratios will vary strongly along the 
line profile, reflecting the changing excitation conditions in the inner region of NGC 1068. 

An important consequence of the moderate 12 CO J=l-0 and the correspondingly 

small (< 0.1) 13 CO J=l-0 optical depths is that the latter transition can be slightly 

super-thermally excited (Goldsmith 1972; de Jong et al. 1975; Leung & Liszt 1976). Indeed 

for the best LVG solution it is T cx (l — 0) > T kin for 13 CO. In this case the sensitivity 

of Rio to the excitation conditions with respect to the R 2 i and R 32 ratios becomes more 

(13) 

pronounced while r 21 ~ 0.8 — 1.0, and r 21 <^ 1.3 can still be fitted with a single gas 
phase of small/moderate 12 CO J=l-0 optical depth (R 10 ^ 10). Later we will see that this 
property of the 13 CO J=l-0 transition remains as one the characteristics of the diffuse gas 
phase in a two-component gas model. 
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5.1. The heterogeneous gas in NGC 1068 

The inability of a single gas phase to account for the line ratios of the 12 CO, 13 CO 
and C 18 isotopes in the starburst region of NGC 1068 suggests at least two distinct 
phases. One phase dominating the 12 CO while the other is responsible for most of the C 18 
emission, and with the 13 CO having contributions from both. 

The notion of a two-component or generally multi-component molecular gas phase 
seems trivial in the face of numerous studies of the Milky Way and other galaxies that reveal 
a large range of temperatures T kin (~ 10 — 100K) and densities n(H 2 )(~ 10 2 — 10 s cm~ 3 ) 
for the molecular gas. This notion is indeed not particularly useful unless one includes 
the spatial scale involved. Examination of the physical conditions of the molecular gas on 
scales of L ^ 1 — 50 pc is bound to yield a very different picture depending on whether one 
looks at Orion-type molecular clouds that harbor ongoing star formation or more quiescent 
ones. Averaging over much larger (L ^ 0.5 kpc) scales "smoothes" out the local density and 
temperature irregularities to a degree that a narrow range of temperatures and densities, 
deduced from molecular line spectroscopy, can fit the observed line ratios. 

It is when such "averaging" fails to converge towards a single narrow set of physical 
conditions that the notion of more than one molecular gas components becomes meaningful 
and necessary. 

Numerous studies of molecular gas in other galaxies (e.g., Knapp et al. 1980; Wall & 
Jaffe 1990; Eckart et al. 1991; Braine & Combes 1992; Aalto et al. 1995) have demonstrated 
a broad differentiation of molecular gas properties with warm (T kin ^ 20 K) and dense gas 
(n(H 2 ) ^ 10 4 cm~ 3 ) lying in the inner ~ 0.5 — 1 kpc, while colder (T kin ^ 10 K), less dense 
(n(H 2 ) ^ 10 3 cm~ 3 ) gas is located further out in the disk. A similar differentiation seems 
to be in place in the Milky Way (e.g., Bally et al. 1988; Stark et al. 1989, 1991; Binney 
et al. 1991; Spergel & Blitz 1992). The properties of the molecular clouds in galactic 
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disks seem to be easily described by a single set of molecular gas properties or equivalently 
the "averaging" of those properties over scales of L ^ 0.5 kpc converges to single type of 
"average" cloud. 

This does not seem to be the case for the molecular gas lying in the central galactic 
regions. Several studies (e.g., Wall & Jaffe 1990; Eckart et al. 1990; Wall et al. 1993; Aalto 
et al. 1995) indicate that even when averaged over ~ 0.5 — 1 kpc, the line emission of the 
various trace molecules does not point towards a single gas phase but at least two. For 
starburst centers this "differentiation" of the molecular gas seems to become more acute 
(Aalto et al. 1995) and is thought to be driven by a variety of mechanisms, namely an 
intense UV field, large turbulent linewidths, tidal forces and cloud-cloud collisions. In this 
picture one phase consists of warm and diffuse and low-density molecular gas originating 
from the heated and disrupted outer envelopes of molecular clouds while the other one 
consists of the more "protected" interiors of such clouds where the gas is denser and colder. 



The wealth of molecular line data gathered for the central starburst region of NGC 1068 
allows for another rigorous test of the aforementioned ideas and a better determination of 
the properties of the two gas phases. 

We denote as (B) a dense, concentrated gas phase where 13 CO and even C 18 may 
have substantial optical depths and as (A) a diffuse phase where t w ~ 1 — 2 for 12 CO 
J=l 0. Then all the observed line ratios can be expressed as a weighted mean of the ratios 
of these two phases. Hence the Rio 8 " 1 line ratio can be expressed as follows: 



5.2. A two-phase model for NGC 1068 





(A) + [f c pi 3 ]Rlo j (B) 

1 + [fcPis] 



(3) 
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where f c is the filling factor of phase (B) relative to phase (A) and p i3 = T R (B)/T R (A) is 
the ratio of the 13 CO J=l-0 brightness temperature for the two gas phases. 

In phase (A) both 13 CO and C 18 isotopes are optically thin (so there is little radiative 
trapping), and since they have the same collisional excitation coefficients and are influenced 
by the same background radiation field it is R^f^A) = [C 18 0/ 13 CO]. Assuming the highest 



such abundance, measured towards the central regions of the Milky Way (e.g., [Wannier 
1980| ; |Hcnkel fc Mauersberger 1993| ; [Langer fc Penzias 1993| ) , means 

RSo 8) (A) = 0.15. Hence, 



R 10 can be expressed as 

R 10 (A) + [f c p 13 ]Rio(B) 
Rio = 7- — f7 ] • (4) 

The rest of the line ratios rj +1 j and rjj+j j (J=l, 2) for 12 CO and 13 CO respectively 
can be expressed as follows: 



_ rj+i j(A) + [f c pi 2 ]rj+i j(B) 

r J+l J - , . r f 1 W 

1 + [tcPial 



,.(18) _ r? + 3 } j(A) + [f^jrgj j(B) 

1 + [fcPis] 



r j+i j — \ mrz i ' \°) 



where p i2 = T R (B)/T R (A) denotes the brightness temperature ratio of 12 CO J=l-0 for the 
two gas phases. The two factors pis and pyi are related by the simple relation 



Rio(B) 



Pis- (7) 



Rio(A) 

We find that the range of the observed line ratios can be reproduced only when 
fcPi3 ^ 1 and f c pi2 1. The latter yields r 2 i ~ r 2 i(A) and r 32 ~ r 32 (A) (Equation 5), i.e., 
these ratios are dominated by the diffuse phase (A). The C 18 J=l-0 optical depth for 
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phase (B) is found to be r} ~ 1 and hence 13 CO and 12 CO in this phase are optically 
thick. Radiative trapping is then enhanced and can thermalize the J < 3, transitions of 
12 CO at densities of n(H 2 ) ~ 10 2 cm -3 , and 13 CO at n(H 2 ) ~ 10 3 cm -3 . non-LTE effects are 
expected to be important for the 13 CO, C 18 transitions arising in the diffuse phase (A). 

In this simple model we assume LTE excitation of the 12 CO, 13 CO J=l-0, 2-1, 3-2 and 
C 18 J=l-0 transitions for the more spatially concentrated gas phase (B), while making 
no such assumption for the "envelope", diffuse gas phase (A). The deduced line ratios for 
the latter are analyzed using LVG modeling with the same parameter space and abundance 
ratios as in our previous single-phase LVG analysis. A good solution for phase (A) is found 
only for r 21 = 0.97, r 32 = 0.74, corresponding to a -30% offset in the BIMA 12 CO J=l 
intensity. The results for both gas phases are summarized in Table 3. 

EDITOR: PLACE TABLE 3 HERE 

5.3. The diffuse gas phase 

From Table 3 we see that this gas phase is characterized by small/moderate optical 
12 CO J=l-0 depths of t 10 <; 1, relatively low density n(H 2 ) ~ 10 3 cm~ 3 and probably 
higher kinetic temperatures than the denser and optically thick phase (B). The 13 CO J=l-0 
transition is now significantly super-thermal, which is responsible for the high Ri = 17 and 
moderate r 21 ~ 1 and ~ 1.4 ratios. This non-LTE effect can occur for a wide range of 
physical conditions as long as T kin ^ 20 K, n(H 2 ) ^ 2 x 10 3 cm -3 and <^ 0.1 (Leung 
& Liszt 1976). Moreover, the higher the T^n the larger the density regime over which 
superthermality occurs, e.g., for Tkm = 60 K and A ~ 10~ 5 (km s _1 pc -1 ) -1 the density 
range is n(H 2 ) ~ 10 3 — 10 5 cm -3 . This includes densities well past the thermalization 
density of the J=l-0 transition (~ 2 x 10 3 cm -3 ). 
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These conditions are more likely to be found in starburst environments rather than 
quiescent ones. There the intense UV radiation and the turbulent motions present will heat 
and disrupt the molecular clouds and thus create the warm diffuse gas phase where 13 CO 
J=l— may become super-thermally excited. The line ratios are then sensitive to (n(H 2 ), A) 
but they generally satisfy r 2 i ~ 0.8 — 1.1, r^ ~ 1.2 — 2.2 and Rio ^ 10, which contains the 
range of the values found for these ratios by Aalto et al. (1995) for starbursts. 

The high value of the R 32 ratio for phase (A) (Table 3) significantly restricts the LVG 
parameter space. This stems from the fact that the conditions for superthermality are 
bounded by the conditions needed to keep a low r 32 (and hence a high R 32 ) which call for a 
low excitation of the J=3 level, i.e., n(H 2 ) < 4 x 10 4 cm -3 and/or T kin < 30 K. However, in 
a starburst environment a wide range of kinetic temperatures T^mfe 30) K can be expected, 
especially for the "exposed" molecular gas in the outer areas of molecular clouds (Aalto et 
al. 1995). Hence the aforementioned restrictions to the excitation conditions of the diffuse 
gas phase become mainly restrictions on its density, namely 10 3 cm~ 3 n(H 2 ) ^ 10 4 cm -3 . 



5.3.1. Gas kinematics and the range of dV/dr 

The other important characteristic of the diffuse gas phase is that is probably not 
virialized. This can be viewed as a direct result of the various mechanisms that create 
this phase by disrupting initially denser and more compact molecular clouds. After the 
"envelope" gas phase is created the energetic environment in which it resides may not allow 
it to settle to self-gravitating structures. Its prevailing kinematic state is important since 
the velocity-integrated luminosity of 12 CO J=l-0 is routinely used as an indicator of gas 
mass in galaxies with one of the key assumptions (among others) being that molecular 
gas consists of an ensemble of virialized clouds, hence linewidth reflects mass. For 12 CO 
emission dominated by non-virialized diffuse gas this will clearly overestimate the gas mass. 
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The LVG approximation is valid when the linewidth Av of the average cloud is 
larger than the thermal one, i.e., Av 3> Av th . For temperatures of T kin ~ 10 — 50 K the 
thermal linewidth is Av t h = 0.1 — 0.2 km s _1 . Hence LVG approximation is expected to 
be valid for linewidths of the order of Av ^ 10 x Av t h = 1 — 2 km s" 1 , which is what 
is usually observed in individual molecular clouds. An upper limit to the linewidth of 
such clouds is set by the fact that their individual linewidths cannot be larger than the 
observed ones in extragalactic systems, which consist of many molecular clouds. For a 
typical extragalactic linewidth of Av G = 300 km s _1 consisting of the linewidths of at 
least ~ 10 clouds, the maximum linewidth is Av max = 0.1 x Av G = 30 km s -1 . This is a 
rather generous upper limit since the velocity dispersion seen in large individual Galactic 
GMCs is of the order of Av = 3 - 8 km s^ 1 (Scoville & Good 1989). For a cloud size of 
L ~ 0.5 — 50 pc, the extreme velocity range of ~ 2 — 30 km s~ x yields a maximum range of 
dV/dr ~ (0.04 — 60) km s _1 pc -1 for the average velocity gradient per cloud. 

If the "average" cloud is virialized then, for a given average gas density, the expected 
range for dV/dr is significantly narrower. For a spherical cloud of radius R the virial 
theorem gives: 



where \i is the mean mass per particle and (n) is the mean number density of the 
cloud. The factor a ~ 0.5 — 2.5 depends primarily on the assumed density profile (see 
Bryant & Scoville 1996). Expressing the last result in astrophysical units gives 



For NGC 1068, the best LVG solution for phase (A) gives (Table 3) (n) = 3 x 10 3 cm" 3 , 
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so for the largest reasonable value of a = 2.5 we estimate (dV/dr) VIR ~ 2 km s 1 pc 1 . 
The same model also gives A = 3 x 1(T 6 (km s" 1 pc" 1 )- 1 . For [ 12 CO/H 2 ] = 10~ 4 this 
yields dV/dr = 33 km s" 1 pc" 1 , which is ~ 17 times larger than the virial value. Hence, 
unless the [ 12 CO/H 2 ] abundance is significantly lower in the diffuse phase, such a large 
value of dV/dr suggests that this phase in not virialized. A similar situation is found for 
the starbursting center of IC 342 by Irwin & Avery 1992, which they interpret in terms of 
a significantly lower [ 12 CO/H 2 ] abundance rather than in terms of non- virialized gas. In 
Ultraluminous Infrared Galaxies, such a gas phase seems to be filling the entire intercloud 
medium (Solomon et al. 1997) in the inner few hundred parsecs. 

In all these cases the 12 CO linewidth will have other significant contributions besides 
the gravitational potential of the gas itself, namely pressure, stellar and other CO-dark 
mass (Downes et al. 1993) and the general potential of the galaxy (Solomon et al. 1997). 



A lower limit for the gas density of phase (B) can be deduced from the information in 
Table 3 and the fact that LTE provides a good description of the excitation of the gas in 
this phase. The optical depth of the C 18 J=l-0 transition can be expressed as follows 



where Tk; n is the kinetic temperature of the gas in (B) phase. 

Assuming that the "average" cloud of this gas phase is a virialized structure means 
that dV/dr is given by Equation 9, so for an abundance [ 12 CO/H 2 ] = 10" 4 we get 



5.4. The dense gas phase 




(10) 
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n(H 2 ) = 2.78 x 10~ A a 



) 



2 ri2 C0 




(11) 



I _ e -5.26/T kin 



c 18 o 



A temperature of Tk in = 10 K and an optical depth of r 10 = 0.7 are the 
minimum values characterizing phase (B). These values and an abundance of 



[ 12 CO/C 18 0] = [ 12 CO/ 13 CO] x [ 13 CO/C 18 0] = 40 x 1/0.15 « 267 yield a gas 



density of n(H 2 ) w 10 4 cm -3 (a = 1.5). This is a minimum value since T kin can attain 
significantly larger values (Table 3) which would yield larger values of n(H 2 ), e.g., for 
T kin = 20 K and = 1.1 it is n(H 2 ) 3 x 10 5 cm -3 . Such densities are comparable to 
the average densities found in the cores of molecular clouds like Orion, M 17 and Cepheus 
by Bergin et al. (1996) which also find good agreement between the virial masses and 
the ones deduced from C 18 observations (Goldsmith et al. 1997). Situations where 
even the rare C 18 isotope may have moderate optical depth and is associated with high 
density molecular gas are known. In these cases C 18 shows good kinematic and spatial 
correspondence with much higher dipole moment molecules like HC3N (Bergin et al. 1996) 
and NH 3 (Fuente 1993). Similarly high densities but with higher average temperatures are 
found also towards the Galactic Center region (Guesten 1989). 

The large abundance ratio [ 12 CO/C 18 0] plays a crucial role in making these two 
isotopes suitable for broadly tracing the two gas phases described previously. Indeed, while 
it is frequently said that molecules with high dipole moments (e.g., CS, HCN) are needed 
in order to trace high density gas, to a certain extent, rare CO isotopes can do the same. 
This happens because unlike the high density tracing molecules, 12 CO has an abundance 
high enough for moderate/high optical depth to arise in the general conditions of the ISM. 
Most of the emission of a transition arises in regions where r ^ 1 and the transition is 
thermalized. Thus 12 CO J=l— can trace a diffuse gas phase since, for t 10 ~ 1 — 10, radiative 
trapping will thermalize this transition at n A (H 2 ) = n crit x (1 — e _Tl() )/r 10 ps 10 2 — 10 3 cm~ 3 , 
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where n crit is the critical density of the J=l-0 transition. On the other hand, assuming 
uniform gas excitation, the less abundant isotope C 18 with [ 12 CO/C 18 0] = 250 — 500 
(Wannier 1980), will attain an optical depth of t w ~ 1 in regions with densities of 
n B (H 2 ) ~ (250 - 500) x n A (H 2 ) £ 2 x 10 4 cm" 3 . 

5.5. The changing physical conditions of the gas 

The different angular resolution of the maps of the various CO transitions prevent us 
from using all the observed line ratios to study the molecular gas within the starburst region 
of NGC 1068 at a common, high angular resolution. Nevertheless, the wide range of physical 
conditions can be revealed from the r 32 and Rio ratios alone and simple one-phase LVG 
modeling. Apart from the fact that both these ratios are sensitive to gas excitation they can 
be estimated at the highest possible resolutions, and hence because of less spatial smoothing 
they can be sensitive to a wider range of physical conditions. For Ri ^ 10 ((Rio) = 14) 
the range of r 32 (~ 0.4 - 0.8) corresponds to T kin > 20K, n(H 2 ) = 10 3 - 10 4 cm" 3 and 
A = 10~ 6 - 10" 5 (km s" 1 pc" 1 )" 1 . For r 32 = 0.57 - 1.1 (-30% offset for 12 CO J=l-0) this 
parameter space becomes even wider. In all cases we find the most dense gas inferred for 
the large GMAs in the western part of the central emission. 

Using all the global line ratios and our simple two-phase model confirms the results of 
the one-phase LVG analysis. For the "low- velocity" CO emission (Vlsr ~ 980 — 1050 km s _1 ) 
corresponding to GMAs in the easternmost part of the central region (Heifer & Blitz 
1995) we find negligible contribution of the dense gas phase (B) to the 12 CO and the 
13 CO emission, and the diffuse phase (A) seems to be the dominant component. For the 
"high- velocity" CO emission (V L sr ~ 1200 — 1280 km s _1 ) originating from the GMAs 
located in the westernmost part there is a significant contribution from the dense and 
optically thick gas phase (B). This picture is further supported by a sensitive map of the 
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HCN J=l-0 which traces dense gas made with the IRAM interferometer (Tacconi et al. 
1994) at a resolution of ~ 4" x 3". This map shows that, apart from the Seyfert nucleus 
itself, the brightest HCN emission originates from the western part of the inner kpc-size 
region of NGC 1068 and is most prominent towards the location of the massive GMA at 
(Aa, AS) « (-8", -15") from the center. 

Our study lacks sufficient angular resolution to shed more light to the physical 
conditions of the nuclear gas (L ^ 100 pc), where warm (T kin ^ 70 K) and dense 
(n(H 2 ) ~ 10 5 cm~ 3 ) gas dominates (Tacconi et al. 1994). As we discussed earlier, even the 
r 32 ratio at ~ 16", is dominated by the extended CO emission from the GMAs throughout. 



5.6. Molecular gas mass in NGC 1068 

The most common method used to find molecular gas mass in extragalactic systems 
relies on the so-called standard galactic conversion factor X G which converts the luminosity 
of the 12 CO J=l-0 line to gas mass, namely 

X G = = 5.0 M (km s" 1 pc 2 )" 1 , (12) 

Leo 

with its numerical value adopted from Solomon & Barett (1991). The 12 CO J=l-0 
luminosity L co is estimated from 

Lco = A s f f T b (a,v) dadv = fi c D 2 / Tj, a) (v)dv, (13) 

JAv JA S JAv 

where Tj^(v) is the area-averaged spectrum (see Figure 5) and D the distance to the galaxy. 

The velocity-integrated brightness of 12 CO J=l-0 from the Kaneko et al. (1989), 
Maiolino et al. (1997) and our measurement (scaled by —30%) give an average of 
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(I) = 88 K km s -1 , over an area of ~ V x 1'. Using X G then gives M(H 2 ) « 6 x 10 9 M . 

Under certain assumptions and taking into account only differences in average density 
(n(H 2 )) and 12 CO J=l brightness temperature (T b ) (see Bryant & Scoville f996 for a 
more complete treatment), the conversion factor X can be expressed as follows (e.g Radford 
et al. 1991) 



From our single-phase LVG analysis we find that X/X G w 1-3, which can yield gas mass 
estimates as high as M(H 2 ) ~ 1.8 x 10 10 M . 

Several studies (Dickman et al. 1986; Maloney & Black 1988; Bryant & Scoville 
1996) attempt to quantify the effects of the molecular gas environment on the value of X, 
especially when the 12 CO emission is optically thick. These studies indicate that, when a 
standard galactic conversion factor is used, the presence of a non-virialized gas component 
will result to an overestimate of gas mass, while significant shadowing of molecular clouds 
in spatial and/or velocity space to an underestimate. Except when the angular resolution 
is high enough to resolve individual clouds, cloud-cloud shadowing may be important only 
towards very quiescent cloud environments unlike the ones expected in a central starburst. 
Therefore the sole most important factor affecting X is the presence of a non-virialized gas 
component where molecular cloud linewidths are no longer good mass indicators. 

A different method is to use rare isotopes of CO like 13 CO or C 18 to deduce molecular 
gas mass. The main assumption here is that these isotopes have small/moderate optical 
depths throughout the volume of the emitting gas. Then the total number of the CO 
isotope molecules at the rotational level J is given by the expression 




(14) 



-27- 



Nj = Tffi' 1 L co(J> J - 1), (15) 

hc ci Aj > j_i 

where Aj j_! and are the Einstein coefficient and the escape probability of the 

J — > J — 1 transition (e.g., pobolev 1960| ; Castor 1970) respectively, and L co (J, J — 1) is 
the corresponding luminosity (Equation 13). A good approximation of the total H2 mass is 
then provided by 



M(H 2 ) = X co (Ni + N 2 + N 3 ) /im H2 , (16) 

where X C o is the abundance of H 2 relative to the particular CO isotope and /im H2 is the 
mean mass per H 2 molecule. Hence the expression for the Xis factor that converts the 
J=l-0 luminosity of the CO isotope to molecular gas mass is as follows 
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where r J+1 j = L co (J + 1, J)/L C q(1, 0). After substitution this equation yields 



X IS = 0.078 [r] 



where we assumed [H 2 / 12 CO] = 10 4 and [r] = [ 12 CO/ x C y O] is the abundance of 12 CO 
relative to the isotope used. 

The assumption about the virialization of the average molecular cloud needed to 
deduce X (Equation 16) is no longer necessary. The omission of Nj for J=0 and J > 3 is 
not expected to have a major effect on the estimated M(H 2 ) for the average conditions 
encountered in molecular clouds. This is suggested by a study of Giant Molecular 



-28- 



Cloud (GMC) cores done by Goldsmith et al. (1997), where they find good agreement 
between the virial mass and the mass deduced using this method and observations of the 
C 18 isotope. They also find that the true mass can be at most ~ 2 times larger than 
the one computed from Equation 16. This occurs for cold (T kin ss 10 K) and diffuse 
(n(H 2 ) < 10 3 cm~ 3 ) gas where the term N becomes important, or warm (T kin > 50 K) and 
dense (n(H 2 ) > 10 5 cm~ 3 ) gas where levels above J=3 become significantly populated. 

In the optically thin case, ~ 1, and X IS is independent of the details of the gas 

excitation conditions and velocity field that determine the value of /Sj+^j's. Of course in the 
optically thick case these become important and more careful analysis is needed to estimate 
the escape probabilities. 

Applying this method for 12 CO or 13 CO transitions, we find a gas mass of 
M(H 2 ) ~ 5 x 10 s M from both isotopes. This value is ~ 10 — 35 times smaller than the 
one deduced from the standard galactic conversion factor X G or its density/temperature 
"corrected" value X. Such a large discrepancy could not be due only to the fact that 
Equation 16 gives a lower limit to the molecular gas mass. Summing up the population of 
all the energy levels can result in a gas mass that is at most ~ 2 times larger (Goldsmith et 
al. 1997), which still leaves M(H 2 ) at least ~ 5 times smaller than the one found from the 
standard galactic conversion factor. 

Two effects that can contribute to this difference in the estimated gas mass are: a) the 
presence of a non-virialized gas component and, b) a large optical depth of the particular 
CO isotope used. The first effect will lead to an overestimate of the molecular gas mass 
whenever a factor that converts 12 CO line luminosity to gas mass is used. Our analysis 
suggests that such a non-virialized molecular gas component indeed exists and dominates 
the 12 CO emission in the central region of NGC 1068. 

The second effect is due to the fact that in an optically thick medium significant 
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amounts of gas mass will contribute only a small fraction of the total observed luminosity. 
This is straightforward to understand when the line formation mechanism is local because 
then the observed emission comes from the outer layer of the cloud where r ~ 1 — 3 and 
a lot of gas mass can "hide" in its interior. The /3j+i,j's factors can, to a certain extent, 
correct for this effect but obviously this correction becomes increasingly unreliable the 
larger the optical depth. 

However in the case of large scale systematic motions within each cloud (and among 
the various clouds themselves), the various optically thick components appear centered 
at different velocities along the line profile and hence are "visible". In such a situation 
significant gas mass can still "hide" provided that the more optically thick component has 
a lower temperature, a smaller velocity gradient and overall volume than the more diffuse 
gas component. This is exactly the type of molecular gas differentiation found for the gas 
in NGC 1068 and generally in starburst centers. A combination of a small spatial/velocity 
filling factor and a lower average temperature can reduce the luminosity contribution of the 
optically thick gas phase to <; 0.1% of the luminosity of the diffuse, warmer gas component 
for equal amounts of gas mass. 

The more accurate description of the molecular g two-phase system allows 

a better estimate of the total gas mass as well as an estimate of the relative gas mass 
contained in each phase, namely 



Mr 

^ = 0.13ftp 13 ) 
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1 + 0.42 + 0.26 r 32 

where Tkm and flio are the average kinetic temperature and escape probability of C 18 
J=l-0 in gas phase (B) and r^ and r^ are the line ratios for the optically thin 13 CO 
(/3 J+ i j ~ 1, J=0, 1, 2) in phase (A). 
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Typical values from Table 3 yield M B /M A ^5 — 10 suggesting that most of the 
mass resides in the dense, optically thick gas phase (B). The estimated total gas mass is 
M(H 2 ) = Ma + Mb ~ 3 x 10 9 M . As expected, this estimate yields a mass that is larger 
(by a factor of ~ 6) than the one found from Equation 18 under the assumption that 13 CO 
is optically thin throughout the emitting volume, yet smaller (by a factor of ~ 2 — 6) than 
the ones estimated from the standard galactic conversion factor. 

5.6.1. Estimates o/M(H 2 ) in starburst environments 

Our study of the IR- luminous starburst in NGC 1068 suggests that in the inner few 
kiloparsecs, the highly excited, non-virialized molecular gas component that dominates the 
12 CO emission can be the cause of a systematic overestimate of molecular gas mass when 
a standard galactic conversion factor is used. In cases of other IR-luminous galaxies with 
powerful central starbursts like Mrk 231 and NGC 7469 the gas mass deduced from the 
standard galactic conversion factor is manifestly overestimated since it yields a mass that is 
comparable or larger than the dynamical mass contained within the same volume (Bryant 
& Scoville 1996; Genzel et al. 1995). While this can still be the result of a particular 
geometry of the CO-emitting gas (e.g., for Mrk 231, Bryant & Scoville 1996) it is also likely 
that the diffuse gas component in their starburst regions is in a highly excited state with 
high temperatures and large non-virial linewidths. In this case the velocity-integrated 12 CO 
J =1-0 luminosity will be significantly larger than in more quiescent environments, without 
this reflecting a proportionally larger molecular gas mass. 

In the particular case of Seyferts, since type 2 are more likely to harbor starbursts than 
type 1 or field spirals (Maiolino et al. 1995), it follows that molecular gas mass deduced 
from their 12 CO J=l-0 luminosity will systematically overestimate their gas mass. 
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6. CONCLUSIONS 

We conducted an extensive set of observations that map the distribution of the 
molecular gas in the inner ~ 1' x 1' of the archetypal Seyfert 2/starburst galaxy NGC 1068. 
The fully sampled maps of 12 CO, 13 CO J=2 1, 3-2 transitions, a 12 CO J=l— map and two 
test spectra of the C 18 J=2-l transition allowed us to probe the physical properties of the 
gas in great detail. Our main conclusions are the following: 

1. A single-phase LVG model yields average conditions that correspond to relatively 
dense molecular gas (~ 10 4 cm -3 ), with a temperature of ~ 20 K and optical depth 

t ~ 1 — 2 for the 12 CO J =1-0 transition. These conditions are not uniform and the 12 CO 
r 32 = (3 — 2)/(l — 0) ratio (at resolution of ~ 16") and earlier interferometric measurements 
of the Rio = 12 CO/ 13 CO and R& 8) =C 18 0/ 13 CO J=l-0 ratios reveal a wide range of 
molecular gas properties in the inner ~ 1' x 1' starburst region of NGC 1068. 

2. A single gas phase with moderate 12 CO J=l-0 optical depths cannot reproduce the 
high C 18 0/ 13 CO (J=l-0, 2-1) intensity ratio and its large spatial variations observed in 
interferometric maps of the J=l-0 transition. A simple two-component model is used to 
describe the state of the "average" molecular cloud as consisting of a dense 10 4 cm -3 ), 
more spatially concentrated component where C 18 J=l-0 has a moderate (r ~ 1) optical 
depth, surrounded by a more diffuse (~ 10 3 cm" 3 ), warmer gas phase where the 12 CO 

J =1-0 transition has r ~ 1 — 2. 

3. The dense and more spatially concentrated phase contains most of the gas mass. 
Since the 13 CO in this phase is optically thick, using this isotope to deduce molecular 
gas mass under the assumption that it has small/moderate optical depth throughout the 
CO-emitting volume underestimates the true mass of this component and hence the total 
gas mass. 
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4. The diffuse gas phase dominates the observed 12 CO emission, it is highly excited 
and probably not virialized. This leads to an overestimate of molecular gas mass when the 
luminosity of the 12 CO J=l-0 line and a standard galactic conversion factor are used to 
deduce gas mass in starburst environments. Since type 2 Seyferts harbor starbursts more 
often than type 1 or field spirals, it follows that the use of the standard galactic conversion 
factor will systematically overestimate the global molecular gas mass in Seyfert 2's. 
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Fig. 1.— The (T mb ) Av maps of 12 C0 (top) and 13 C0 (bottom) of the J=2-l transition, 
with Av = 950 — 1300 km s -1 , at a resolution of 0; = 24". For 12 CO the rms noise is 
<7i2 ~ 0.025 K, the contours are (3, 5, 7, 9, 11, 13, 15, 17, 19, 21) x 0i 2 , and the greyscale 
range is 0.075-0.535 K. For 13 CO the rms noise is 013 0.003 K, the contours are (3, 5, 7, 
9, 11, 13, 15)x0 13 , and the greyscale range is 0.009-0.048 K. 

Fig. 2.— The (T mb ) Av maps of 12 CO (top) and 13 CO (bottom) for the J=3-2 transition, 
with Av = 950 - 1300 km s _1 , at a resolution of 6 { = 16". For 12 CO the rms noise is 
012 ~ 0.020 K, the contours are (3, 5, 7, 9, 11, 13, 15, 17, 19, 21) x 012, and the greyscale 
range is 0.060-0.456 K. For 13 CO the rms noise is 013 ss 0.002 K, the contours are (3, 5, 7, 
9, 11, 13) x 013, and the greyscale range is 0.006-0.026 K. 

Fig. 3. — Channel maps of the 12 CO J=3-2 emission (contours) overlayed on the 12 CO 
J=l-0 emission (greyscale) at a common spatial resolution of 16" and velocity resolution of 
A v chan ~ 8 km s -1 . The average rms noise of the 12 CO J=3-2 map is 0t(3 — 2) = 0.065 K, 
the lowest level contour is 3 T (3 — 2) and the contour interval is 2 T (3 — 2). The rms noise 
of the 12 CO J=l-0 map is T (1 - 0) = 0.08 K and the greyscale range is: 0.25-2.50 K. The 
line ratio is estimated as rf^^ = (T m b(3 — 2))r l / (T m b(l — 0))r l for a radius of Rl <^ 30" 
from the map center and it is reported at the bottom left corner of each channel map. The 
beam is shown at the bottom right corner of the first channel map. 

Fig. 4. — Spectra of C 18 0, 13 CO J=2-l in two positions The offsets from the map center are 
the following A: (Aa, AS) = (8",0"), B: ( Aa, AS) = (-8", -3") 

Fig. 5. — Area-averaged spectra of 12 CO, 13 CO for J=3-2, 2-1, 1-0, over a circular area 
with radius R = 30". 
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Table 1 

NGC 1068: The observations 



Run 


Receiver 


Spectral line 


Av c 


#HPBW d 


A# s c 


Map size 




(T S y S ) 


K) b 


(Av chan ) c 


(<5#rms) d 


(T int ) c 


(# of points) f 


1994, Jan. 2-8 


A2 


12 CO J=2-l 


900 


21" 


10" 


100" x 60" 


1994, Nov. 26-29 


(350-450) 


(230.538) 


(0.81) 


(4") 


(1-2) 


(85) 


1994, Jan. 2-8 


A2 


13 CO J=2-l 


941 


21" 


10" 


80" x 60" 


1994, Nov. 26-29 


(400-500) 


(220.398) 


(0.85) 


(4") 


(10-20) 


(57) 


1997, Nov. 27 


A2 


C 18 J=2-l 


945 


22" 








(500) 


(219.560) 


(0.85) 


(4") 


(90, 110) 


(2) 


1994, Nov. 26-29 


B3i 


i2 C J=3-2 


600 


14" 


7" 


80" x 60" 


1996, Jan. 16-28 


(750-1000) 


(345.795) 


(0.54) 


(3") 


(2) 


(113) 


1997, Nov. 27-Dec. 1 


B3 


!3C0 J=3-2 


900 


14" 


7" 


60" x 50" 


1997, Dec. 27-29 


(400-500) 


(330.588) 


(0.57) 


(3") 


(10-15) 


(87) 



a Typical system temperatures for the corresponding receiver in Kelvin. 
b The rest frequency of the spectral line in GHz. 

c The corresponding velocity range Av, and resolution Av c h an in km s _1 . 

d The HPBW of the gaussian beam, and the rms pointing error S9 rms . 

e The sampling interval A^s, and the integration time T int per point in minutes. 

f The number of points corresponding to each map. 
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Table 2 

NGC 1068: The global line ratios 



12 CO(3-2) 
12 CO(1-0) 

( r 3 2 ) 


12 CO(2-l) 
12 CO(1-0) 

M 


12 CO(3-2) 
13 CO(3-2) 

(R 32 ) 


12 CO(2-l) 
13 CO(2-l) 

(R2l) 


12 CO(1-0) 
13 CO(1-0) 

(Rio) 


C 18 O(l-0) 
13 CO(1-0) 

(r!o 8) ) 


C ls O(2-l) 
13 CO(2-l) 

(rS 8) ) 


0.52 ±0.17 


0.68 ±0.23 


14 ±3 


10 ±2 


14 ±2 a 


0.30 ±0.10 b 


0.23 ±0.04 c 



a OVRO measurement (Papadopoulos, Seaquist & Scoville 1996), Young & Sanders (1986). 
b OVRO measurement (Papadopoulos, Seaquist & Scoville 1996). 
c The average of the two values measured with the JCMT. 
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Table 3 

NGC 1068: Two-phase model 



Gas phase (B) a 

(18) r rp 
TlO ) r cPl3) J-kin 


Line ratios (A) b 

f21, r 32, R-10, R-21, R-32 


LVG parameters (A) c 
T kin , n(H 2 ), A 


Excitation (A) d 

F (13) (12) 
^10 ) r 10 


D 5 1 35 1 5 

_L . '.J '.J • X xj 


97 74 32 15 31 
(1.18, 0.73, 27, 23, 31) 


35 3x1 3 1x1 0~ 6 
X 2 = 3.2 


3 8 29 


7 70 1 D 


97 74 23 1 3 20 
(1.05, 0.76, 17, 13, 20) 


35 3x1 3 3x1 -6 

■ J '.J . UA1U j <JA1U 

X 2 = 3.2 


2 8 94 


9 50 15 


97 74 20 12 20 
(1.05, 0.76, 17, 13, 20) 


35 3xl0 3 3xl0~ 6 
X 2 = 0.95 


2.8, 0.94 


1.1, 0.40, 20 


0.97, 0.74, 19, 12, 20 
(1.05, 0.76, 17, 13, 20) 


35, 3xl0 3 , 3xl0~ 6 
X 2 = 0-53 


2.8, 0.94 


1.5, 0.30, 20-40 


0.97, 0.76, 18, 12, 19 
(1.05, 0.76, 17, 13, 20) 


35, 3xl0 3 , 3xl0~ 6 
x 2 = 0.37-0.45 


2.8, 0.94 


2.5, 0.25, 15-40 


0.97, 0.74, 17, 11, 17-18 
(1.05, 0.76, 17, 13, 20) 


35, 3xl0 3 , 3xl0~ 6 
X 2 = 1.2-0.8 


2.8, 0.94 



a The physical conditions of the gas in phase (B): is the optical depth of the C 18 J=l-0 
transition, Tk; n is the kinetic temperature and f c pi3 is defined in Equation (3). 

b The line ratios estimated for gas phase (A), the ones in the parenthesis are the values obtained 

r 1 2 

from the LVG model with the minimum x 2 , where \ 2 = J2i R-(i) ~~ Robs and Rp), R Q b s are 
the model and observed line ratios respectively, with o\ being the corresponding la uncertainty. 

c The LVG parameters (see text) corresponding to the best fit of the line ratios of gas phase (A). 



d E ( 13 ) 



10 



Tcxc/Tkin for 13 CO J=l-0 and t\^> is the optical depth of the 12 CO J=l-0 transition 



.(12) • 



for gas phase (A). 
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Figure 4 



Figure 5 
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